function [p] = prices(q_ij,beta,delta,s_i,cobb_par,E,M,r,K_obs,c_inv,sf)
%UNTITLED Summary of this function goes here
%   Detailed explanation goes here
etilde_ij = entrants(q_ij,s_i,cobb_par,M,delta,sf);
Gunc = etilde_ij./(E - nansum(etilde_ij,2));
e = inverse_m(nansum(q_ij,2),s_i,cobb_par,M,sf);
lambda_e = nansum(q_ij,2)./e;

tt = (1-beta.*delta.*(1-lambda_e))./lambda_e;
p = r - (log(Gunc)+K_obs).*tt;
p = p - c_inv./lambda_e;
if sf==1
    p = r - (log(Gunc)+K_obs) - c_inv;
end
end

